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ABSTRACT 

The analytical phase-space distribution function (DF) of spherical self-consistent 
galaxy (or cluster) models, embedded in a dark matter halo, where both density 
distributions follow the Hernquist profile, with different total masses and core radii 
(hereafter called HH models), is presented. The concentration and the amount of the 
stellar and dark matter distributions are described by four parameters: the mass and 
core radius of the reference component, and two dimensionless parameters describing 
the mass and core radius of the halo component. A variable amount of orbital 
anisotropy is allowed in both components, following the widely used parameterization 
of Osipkov-Merritt. An important case is obtained for a null core radius of the halo, 
corresponding to the presence of a central black hole (BH). 

Before giving the explicit form for the DF, the necessary and sufficient conditions 
that the model parameters must satisfy in order to correspond to a consistent system 
(i.e., a system for which each physically distinct component has a positive DF), are 
analytically derived. In this context it is proved that globally isotropic HH models are 
consistent for any mass ratio and core radii ratio, even in the case of a central BH. 
In this last case the analytical expression for a lower limit of the anisotropy radius 
of the host system as a function of the BH mass is given. These results are then 
compared with those obtained by direct inspection of the DF. In the particular case of 
global isotropy the stability of HH models is proved, and the explicit formula for the 
differential energy distribution is derived. Finally, the stability of radially anisotropic 
HH models is briefly discussed. 

The expression derived for the DF is useful for understanding the relations between 
anisotropy, density shape and external potential well in a consistent stellar system, 
and to produce initial conditions for N-body simulations of two-component galaxies or 
galaxy clusters. 



Subject headings: galaxies: elliptical - stellar dynamics - dark matter 
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1. Introduction 



Recent ground based observations (Miller et al. 1995), and with the Hubble Space Telescope 
show that the spatial luminosity distributions of elliptical galaxies approach a power-law form 
p(r) oc r~ 7 at small radii, with < 7 < 2.5 ( Crane et al. 1995 , Jaffe et al. 1994 , Ferrarese et 



al. 1994 , |Lauer et al. 1995 , [Kormendy et al. 1995 , Byun et al. 1996 ). These findings increase 
considerably the interest of theorists in the study of cuspy models. Two important families 
of spherical dynamical models with a central divergent density that have been explored so far 
are the R 1 /" 1 models and the so-called 7-models. The dynamical properties of models whose 
surface brightness distribution follows the i? 1//m -law, introduced by Sersic (1968) as a natural 
generalization of the de Vaucouleurs law ( |de Vaucouleurs 1948| ), have been extensively studied 
(Ciotti 1991, Ciotti and Lanzoni 1996). Particularly, their deprojected density increases toward 
the center as r -( m - l )/ m f or rn > 1; unfortunately two major problems afflict these models: their 
deprojected density cannot be expressed analytically in terms of known functions, and no galaxies 
with 7 > 1 can be accurately modeled in their central regions. The family of the 7-models, in some 



way anticipated by Hernquist ( Hernquist 1990 , hereafter H90), has been widely explored ( Dehnen 
1993 , Carollo 1993j , Trcmainc et al. 1994 ) and it represents a generalization of the well known 



Hernquist (H90) and Jaffe ( |Jaffe 19*83 ) density distributions. As shown by the previous authors, 
many of the dynamical properties of the 7-models can be expressed analytically. In particular the 
Hernquist model (hereafter H model) in projection well resembles the de Vaucouleurs law, and an 
exhaustive analytical investigation of its properties is possible (H90). 

It is now accepted that a fraction of the mass in galaxies and clusters of galaxies is made of 
a dark component, whose density distribution differs from that of the visible one. The shape of 
the dark matter distribution is not well constrained by observations, but numerical simulations of 
dissipationless collapses seem to favor a peaked profile, consistent with the scale-free nature of the 



gravitational field ( Dubinski and Carlberg 1991 , White 1996 , and references therein). From these 
considerations it follows that the obvious generalization of the one-component spherical models 
(the dynamicists zero-th order approximation of real galaxies) is not only in the direction of the 
actively developed modeling of axisymmetric and triaxial systems [see, e.g., de Zeeuw (1996) for a 
recent review] but also in the study and construction of two-component analytical models, a field 
far less developed. From this point of view the zero-th order approximation of realistic galaxies is 
the construction of analytical spherically symmetric two-component galaxy models. 

When studying a dynamical model (single or multi-component) the fact that the Jeans 
equations have a physically acceptable solution is not a sufficient criterion for the validity of the 
model: the essential requirement to be met by any acceptable dynamical model is the positivity 
of the DF of each physically distinct component. A model satisfying this minimal requirement 
(much weaker than the model stability) is called a consistent model. Two general strategies can 
be used to construct a consistent model or check whether a proposed model is consistent: the "/ 



to p" and the il p to /" approaches ( Binney and Tremaine 1987 , Chap. 4, hereafter BT87). An 
example of the first approach is the extensive survey of two-component, spherical, self-consistent 
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galaxy models carried out by Bertin and co-workers (Bertin et al. 1992). They assume for the 
stellar and dark matter components two distribution functions of the form (and so positive by 
choice) (Bertin and Stiavelli 1984). The main problem with this approach is that generally the 



spatial density is not expressible in terms of known functions, and so only numerical investigations 
are feasible. 

In the second approach the density distribution is given, and assumptions on the model 
internal dynamic are made, making the comparison with the data simpler. But the difficulties 
inherent in the operation of recovering the DF in many cases prevent a simple consistency analysis. 
In particular, in order to recover the DF of spherical models with anisotropy two techniques 
have been developed from the original Eddington (1916) method for isotropic systems: the 
Osipkov-Merritt technique ( psipkov 197E , Merritt 1985| , hereafter OM), and the case discussed by 



Cuddeford and Louis ( Cuddeford and Louis 1995| , and references therein). Examples of numerical 



application of the OM inversion to two-component spherical galaxies can be found in the literature 
(see, e.g., |Ciotti and Pellegrini 1992] , hereafter CP92; |Carollo et al. 19"95 ). For axisymmetric 



systems recently a new inversion technique, less restrictive than the classical ones (Lynden-Bell 



1962, |Hunter 1975| , Dejonghc 1986j ), has been found ( Hunter and Quian 1993] ). If one is just 



interested in the consistency of a stellar system the previous methods give "too much", i.e., give 
the DF. A simpler approach, at least for spherically symmetric multicomponent systems with OM 
anisotropy - as the case discussed in this paper - is given by a method described by CP92, that 
requires information only on the radial density profiles of each component. 

Despite all these efforts, a small number of one-component systems in which both the spatial 
density and the DF are analytically known is at our disposition, and in the more interesting case 
of two-component systems only the very remarkable axisymmetric Binney-Evans model is known 
(Binney 1981, Evans 1993j ). It is therefore of particular interest the result here proved that also 



the DF of HH models with OM anisotropy is completely expressible in an analytical way. This 
family of models is made by the superposition of a stellar and a dark matter distribution both 
following the Hernquist profile, with different total masses and core radii. The concentration and 
the amount of the stellar and dark matter distributions are described by four free parameters, 
and the orbital anisotropy is allowed in both components, following the OM prescription. A 
particularly interesting case is obtained for a null core radius of the " halo" , so mimicking a central 
BH. The study of HH models is also useful for many different reasons: to provide an analytical DF 
for a two-component cuspy system for which the analytical solution of the Jeans equations is also 
available ( Ciotti et al. 1996| ); to investigate the role of anisotropy and mass distribution of each 



component in determining the positivity of their DF; to compute in an accurate and "easy" way 
the model line-profiles, to arrange initial conditions for numerical simulations of two-component 
systems. 

In Section 2 I briefly review the method presented in CP92, formulating it in a way suitable 
for its application to the present problem. Then in Section 3 I introduce the HH models, and 
use the previous method to discuss the limits imposed on their parameters by the positivity of 
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the DF of the two components. It is proved that globally isotropic HH models are consistent 
for any mass ratio and core radii ratio, even in the case of a central BH. In this last case the 
analytical expression for a lower limit of the anisotropy radius of the host system as a function of 
the BH mass is given. In Section 4 I derive the DF for HH models, and their differential energy 
distribution in the case of global isotropy; some velocity sections of the DF are also shown. In 
the case of a dominant halo it is found that the DF of an HH component can be expressed only 
through elementary functions. A particular case - corresponding to a Hernquist model with a 
central BH - is extensively discussed. In Section 5 the exact boundary of the region of consistency 
in the parameter space is obtained using the DF, and the results are compared with those given 
in Section 3. In the same section the stability of globally isotropic HH models is proved, and a 
discussion on the stability of the anisotropic case is given. Finally in Section 6 the main results 
are summarized. 



2. The Consistency of Multi-Component Systems 

More important than the construction of its spatial and projected velocity dispersion profile 
is checking whether a galaxy (or a galaxy cluster) model is consistent, i.e., is described by a 
DF everywhere non-negative. If a system is described as a sum of different density components 
Pk, then each f k must be non negative. This requirement leads us to introduce the concept of 
consistent multi-component decomposition of a system, as discussed in CP92, together with the 
main theorem used here. This theorem permits us to check whether the DF of a multi-component 
spherical system where the orbital anisotropy of each component is described by the OM anisotropy 
is positive, without calculating it effectively. In the OM formulation the radially anisotropic case 
is obtained as a consequence of assuming / = f(Q) with: 

«-«-^. <» 

where £ and L are respectively the relative energy and the angular momentum modulus per unit 
mass, f(Q) = for Q < 0, and r a is the so-called anisotropy radius. With this assumption the 
models are characterized by radial anisotropy increasing with the galactic radius, and in the limits 
r a — > oo the velocity dispersion tensor is globally isotropic. For a multi-component spherical 
system, the simple relation between energy and angular momentum prescribed by equation (1) 
allows to express the DF of the k-th component as: 

, ((D , J d_ rQ* dg k d^ T _ _j_ rQ* d 2 6k d^ T 

MW V8tt 2 dQ k Jo dV T VQk - V8tt 2 Jo d^\ VQk- *t' 1 ' 



where 

Qk = Pk x I 1 + — ) , (3) 



* T ( r ) = £* k (r) is the relative total potential, Q k = £ - L 2 /2rl k , and < Q k < * T (0). The 
second equivalence in equation (2) holds for untruncated systems with a finite total mass (see, 
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BT87, p. 240), as the HH models here discussed. The original theorem, given in CP92 in terms 
of the model radius, is here formulated in terms of the relative potential ^ k of the investigated 
component, since this formulation makes the treatment of HH models easier. 



Theorem : Necessary condition for the non negativity of f k given in equation (2) is: 

dg k 



dV k 



> 0, < tf k < tf k (0). 



(4) 



If this necessary condition is satisfied, a strong sufficient condition (SSC) for the non negativity of 
fk is: 



d 
d^k~ 



/d^T^- 1 
\d^J 



dgk ( d^T 
d¥i 



> 0, < * k < tf k (0). 



(5) 



Proof: See CP92. 

Note that a weak sufficient condition (WSC) 



d 



dg k ( d$r 
dV k \dV k 



> 



(6) 



is obtained in a much easier way using the last expression in equation (2), requiring that > 0. 
The WSC is obviously better suited than the SSC for analytical investigations, due to the absence 
of the weighting square root of the total potential. 

Some considerations follow looking at the previous conditions. The first is that the violation 
of the necessary condition [eq. (4)] is connected only with the radial behavior of p k and the value 
of r ak , and so this condition is valid independently of any other interacting component added to 
the model. Even when the necessary condition is satisfied, f k can be negative, due to the radial 
behavior of the integrand in equation (2), which depends on the total potential, on the particular 
p k , an d on r ak : some permitted values of r ak satisfying the necessary condition must be discarded. 
Naturally, the true critical anisotropy radius is always larger than or equal to that given by the 
necessary condition, and smaller than or equal to that given by SSC and WSC. The previous 
analysis has been performed obtaining analytical or numerical limits on r a for some widely used 
models: for example, in CP92 the King (1972), de Vaucouleurs (1948), and quasi-isothermal 
density distributions were discussed, and more recently the Jaffe (1983), Hernquist (H90) and 
Plummer (1911) distributions ( Ciotti et al. 1996| ). 



A complete analysis of the HH models is given in the next paragraph, deriving the analytical 
constraints on 7" a for 1) an H model, 2) an H model with a central BH, and finally 3) showing the 
consistency of the general HH model in the case of global isotropy. 
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3. The HH Models 

In the study of the HH models, some simplification arises from the fact that both components 
are described by the same functional form. So, I do not distinguish between "stars" and "dark 
matter" , but simply between a reference system and a halo system. The mass M and the core 
radius r c of the reference system are the normalization constants, so that its density distribution 
is: 

p(r) = PN p(s) = s(1 ^ s)3 , (7) 

where s = r/r c and p^ = M/2-rrr^ (H90). The halo density is described by another Hernquist 
distribution, of mass Mh = pM and core radius r^ = (3r c : 

, s. P^p(3 fQ , 

Ph(r) = WTW (8) 

Then the HH density profiles are fully determined by fixing the four independent parameters 
(M, r c ,p, (3), with < p and < (3. Note that for \i = the HH models reduce to the H model, 
and that for (3 < 1 the halo density is more concentrated than the reference component. 

The fundamental ingredient in recovering the DF is the relative potential that for the 
reference component is 

*N 

*(r) = * N tf (s) = t-j-, < * < 1, (9) 

and for the halo 

*h(r) = ^N^h(s) = (10) 

with * N = GM/r c (H90). Note how for (3 = the halo potential is that of a BH of mass M h 
placed at the center of the reference system. As it will become clear in §4, a fundamental property 
of HH models is that their total potential can be expressed as a simple function of the reference 
potential, namely 

§ T = § x ( 1 + - , V 6 = /3 — 1 > -1, 1 (11) 



1 + 6* 

where the interval < * < 1 is monotonically mapped onto < *t < 1 + p/ (3. 



3.1. The Necessary and Sufficient Conditions for the H Model 

Here I apply first the necessary condition to the H model to determine the critical anisotropy 
radius such that a higher degree of radial OM anisotropy produces a negative DF, no matter what 



lr The case 6 = (i.e. (3 = 1) is discarded from now on, it being the case of an H model of total mass M(l + /i), 
fully discussed in H90. 
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kind of halo density distribution is added. The first step is to express the modified density g [see 
eq. (3)] as a function of the relative potential obtaining: 



i - ^ 



i + 



(12) 



where s a = r a /r c is the dimensionless anisotropy radius. 



As shown in Appendix A [eqs. (Al)-(A2)], the necessary condition can be treated analytically, 
requiring 

(3§ M - 2) (1 - * M ) 



\ (4 - 3* M ) 



~ 0.128, 



M 



(13) 



where \I/m is the value of the potential for which the r.h.s. of equation (Al) is maximum. In Fig. 
1 the solid line at the bottom represents the derived lower bound for the anisotropy radius. 

Moving to the study of the sufficient conditions, I apply first the SSC. The calculations can 
be performed analytically, and the result is: 



s a > 



3(5§m~2)(1-§m) 3 1 
15^-39^ + 28 * M 



~ 0.25, 



(14) 



[eqs. (A3)-(A5)], a limit obviously higher than that obtained from the necessary condition. The 
true limit on s a for the H model is within the two previous values, and in fact its value determined 
directly from the DF is ~ 0.202 (see §4.2). In Fig. 1 this value is represented by the solid line 
in the middle. The result of the application of the WSC to the H model will be obtained as a 
limiting case of the more general analysis done in the next paragraph. 



3.2. Sufficient Conditions for the H+BH and Isotropic HH Models 

In order to proceed further with this analytical discussion, and to consider the more 
complicated case of the presence of the halo, we need to use the WSC rather than the SSC, due 
to the presence of the square root of the total potential in the SSC. I now prove two important 
results: 

1. in the case of a H+BH model (i.e., an H model with a central BH of mass fiM) the WSC 
permits us to recover analytically a minimum value of the anisotropy radius as a function of 
the central BH mass; 

2. globally isotropic HH models can be consistently constructed for any value of 
Particularly this means that a globally isotropic H model can consistently host a BH of any 
mass at its center. 
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In the case of a central BH the WSC prescribes that: 



3^mQ») 4 ~ 10^m(^) 3 + 6(2 + ^m(v) 2 - 6(1 + ^mG») + (1 + m) 1 
\ 3* M (/x) 2 - 8^ m (//) + 6 + 6// * m (m) 

[eqs. (A6)-(A10)], and corresponds to the uppermost solid line in Fig. 1. As intuitive, an increase 
of the BH mass produces an increase of the critical anisotropy radius, i.e., a decrease in the 
maximum value allowed for the radial anisotropy. Two interesting consequences can be obtained 
performing the asymptotic analysis of the previous equation, for jjl — > oo and jjl — > 0. 

The asymptotic behavior of equation (15) for \i — > oo is s a (/x) > l/\/2 + 0(1/ fi), and so just a 
slight reduction in anisotropy with respect to the one-component model is required by the presence 
of a BH of any mass; and the finite value of s a (oo) implies that a globally isotropic model can host 
a central BH of any mass. The limiting case of the previous analysis for n = is the WSC applied 
to the H model, and completes the discussion given in §3.1. The limit for /j, —>■ of equation (15) is 



(3*M-1)(1-*M)3J_ 

\ 3«& - 8* M + 6 * M 

and the comparison of this value with the limit on the anisotropy radius derived from the SSC 
(sa > 0.25) is instructive: the difference is due to the weight factor of the square root of the 
potential, contained in the SSC and absent in the WSC. The numerical evaluation of the SSC in 
the case of a BH is very easy, and the result is plotted in Fig. 1 (dot-dashed line). 

A second interesting case for which the WSC can be treated analytically is that of a 
completely isotropic component of the HH model: in Appendix A it is shown that its DF is 
positive for any choice of (/J.,/3). Incidentally, this gives another proof that it is always possible to 
couple consistently a BH of any mass with a globally isotropic H model, in accordance with the 
previous result. 



4. The DF of HH Models 

After the preliminary discussion we can now proceed to the explicit recovering of the DF. Due 
to the fact that both density components of the HH model are described by the same functional 
form, it suffices to compute the DF for the reference component with generic (fi,f3), and then also 
the DF for the halo component is easily recovered. As for the density and the potential, also for / 
it is useful to work with dimensionless functions. So, for the reference component / = /n/(Mj P', Q) 
with /n = / on^'n 3//2 ' an d < Q = Q/^n < 1 + fi/P- For the halo component a similar functional 
form holds, where the two dimensional constants are now the mass and the core radius of the 
halo, and the two dimensionless parameters are = M/M^ = 1/fi and (3h = r c /rh = 1/(3. The 
DF for the halo component is then obtained from that of the reference component changing the 
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normalization constant /n to fab., substituting the dimensionless parameters (fi,(3) with their 
inverses, and re-scaling the parameter Q to the halo central potential. 

The easiest way to compute the DF is to use the first of the identities in equation (2). For the 
evaluation of the integral one would be tempted to express q(^?t) eliminating the radial coordinate 
from the modified density and the total potential: this can be formally done, but the resulting 
expression for the radius involves a quadratic irrationality, that after insertion in equations (3) 
and (7) produces an intractable expression. Here I follow another approach: instead of eliminating 
the radius, the variable of integration is changed from the total potential to the potential of the 
reference component. This is equivalent to a remapping of the domain of definition of the DF, 
from the range of variation of to the range of variation of and leads to introduce a new 
parameter q, defined from equation (11) as: 



Q = qX { l + Ttb q ^) , °- Q - 1 - (l7) 

With this change of variable, the DF is given by equations (B1)-(B3), but as shown subsequently 
in Appendix B, this is again not the natural parameterization for /, that is finally obtained 
defining the variable: 

I 2 = 1 + bq. (18) 

After normalizing to the dimensional scales of the reference component, its DF can be formally 
written as: 

-l 



HQ) = HQ) + 



UQ) h dQ\ d 



^(0 + ^ 



i = r L (Q), (i9) 



»a v^vr 2 \dl J all 
where the subscripts refer to the isotropic and anisotropic parts of the DF respectively, and 

d A - 2 1 -^ (20) 

The sign ± in equation (19) corresponds to the case b > and — 1 < b < respectively, and F^(l) 
and F^(l) are given in Appendix B. Note how with the followed procedure f(Q) results from the 
elimination of the parameter q between equations (17) and (18)-(19). A first (but algebraically 
cumbersome) check of the derived formula is obtained evaluating analytically its limit for n = 0, 
and recovering the DF of the H model given by H90. The lengthy proof is not given here. In the 
general case of /j, > 0, a check is obtained by confrontation of the analytical DF with that derived 
by direct numerical inversion of equation (2): the two families of curves are indistinguishable, with 
percentual errors smaller than 10~ 5 . 

In Fig. 2 (upper panel) some DFs are shown in the case of global isotropy. Their main 
characteristic is that when the halo is more extended than the reference component (dotted lines) 
they are more peaked than the DF of the H model (solid line). On the contrary, for haloes more 
concentrated than the reference component (short-dashed lines), the DFs are flatter. A particular 
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case is that corresponding to a central BH (long-dashed lines): the DFs are positive - as discussed 
in §3.2 - but not monotonically increasing near the model center. In the lower panel of Fig. 2 the 
DF for the same models above are shown, but in this case the anisotropy radius is fixed to s a = 1. 
The same qualitative comments as in the isotropic case apply. 

An important feature of Fig. 2 is the radically different behavior of / for q — ► 1: while for 
any finite core radius of the halo the DF of the reference component diverges near the center, in 
the presence of a central BH the DF converges, i.e., the DF for the BH case cannot be obtained 
directly as a limit for (3 — > of the DF with (3 > 0. The discontinuity in the DF behavior at high 
energies is due to the coefficient of the function X2 in equation (B12): the limit for (3 —* of the 
product between the coefficient and X2 does not vanish, on the contrary, for (3 = this term is 
zero. The reason for this behavior is that the halo potential is uniformly continuous on s E [0, 00 [ 
for > 0, but for (3 = the uniform continuity is lost, and so the equivalence of the limit for 
(3 —* before and after integration is not guaranteed anymore. 

A more direct way to look at this discontinuity is to perform an asymptotic expansion of the 
DF for q — ► 1. A brute force (and highly error prone) approach would be the expansion of the 
expressions given in Appendix B for q — > 1, but a cleaner asymptotic expansion can be obtained 
instead expanding directly equations (B1)-(B3): for (3 > the leading term is 

~ 3 [! + (/?- l)qf j3 

JKA) SV2W+[l + {P-l)q] 2 ^+W{l-qf/^ q 1 j 

The order of divergency, also in presence of the halo, is equal to that given in H90, and assuming 
/i = or (3 = 1 the correct asymptotic formula for the H model is recovered .0 

Consistently with the previous discussion the asymptotic formula for q — > 1 in the case of a 
central BH cannot be obtained as the limiting case of the previous formula when (3 = 0, due to a 
change in the order of the singularity in equation (B2). The correct treatment in this case gives: 

/(q)~ 2 * 9-1, (22) 



and the convergence to of the DF is proved. Finally, note how the presence of anisotropy does 
not affect the behavior of the DF for high relative energies, as dictated by equation (12), where it 
can be easily seen that the divergence of the modified density near the model center - the leading 
term in the asymptotic expansion of the DF - is only due to its isotropic part. The asymptotic 
formulae (21)-(22) have been checked numerically for many choices of (/i, /?) by direct comparison 
with the DF, and the agreement is excellent. 



All the higher-order terms in equations (21)-(22) can written explicitly, but their expressions become increasingly 
complicated. Note that equation (21) of H90 is incorrect, as can be shown expanding equation (17) there. Note also 
that the parameter q used in H90 is the square root of the parameter used in this paper. 
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4.1. Velocity Sections of the DF 

The representation of the DF as function of the integral of motion Q is not easily interpreted 
when orbital anisotropy is present. More intuitive are the velocity sections of /, i.e., for some fixed 
r the distributions of v r and v t , the radial and tangential velocity components, respectively. In 
Fig. 3 the sections f(r,v T ,0) and /(r, 0, v t ) for the H model with s a = 1 are shown. Obviously, 
for global isotropy the two distributions are equal, and for anisotropic systems they become more 
and more similar moving from radii greater than r a to smaller radii, according to the radial trend 
of anisotropy implied by the OM parameterization. As expected, for a finite value of r a , moving 
outward the tangential orbits become more and more de-populated, and this is compensated by 
an increase in the number of high velocity radial orbits. The presence of a massive diffuse halo 
does not alter qualitatively the distributions, that maintain the same aspect. The only peculiar 
characteristic appears again in the case of the central BH, when the maximum of the distribution 
is placed off-center: this fact reflects the central behavior of the DF described previously. As a 
final and general remark on these sections, one can note how they depart appreciably from a 
Maxwellian distribution. 

4.2. The Differential Energy Distribution of HH Models 

In the case of global isotropy the differential energy distribution dM/d£ for each component 
of the HH models can be derived analytically as a function of the parameter I. As shown by BT87 
(p.243), dM/d£ = f{£)g{£), where 

g(£) = lev^TrVG^/) (23) 

is the density of states, and = r^J/ 2 . The explicit form for G ,± (/) is given in Appendix C, 
and the sign ± correspond to the case b > and — 1 < b < respectively. Again, a first check of 
the derived formula can be obtained evaluating analytically its limit for [i = 0, and recovering the 
g{£) as given by H90; as in the case of the DF, the proof is cumbersome, and not shown here. 
When also a halo is present, the check is obtained by direct comparison with the density of states 
derived by numerical integration of equation (C2). Over the whole energy range the curves are 
indistinguishable. dM/d£ is plotted in Fig. 4: when the halo is more extended than the reference 
component (j3 > 1, dotted lines), dM/d£ is not steadily increasing in the outer part of the system 
anymore. 

4.3. The DF for Halo Dominated HH Models 

In the case of a halo dominated model, i.e., when the self-gravity of the reference component 
is negligible, the DF and dM/d£ can be expanded for ^ — > oo, and the resulting expressions are 
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combinations of elementary functions of I. For brevity, and due to its major importance, the 
explicit form of the functions entering the DF are given only. 



4.3.1. The Case (3 > 1 
In this case, and are given by equations (B8)-(B9) with: 

C 6 (0 = Aa r ccos(I)+v^(^ + ^ + i). (24) 
C 4 (l) = § arccos (}) + s/P~l ( J5 + Jj) , (25) 

1 /1\ \/l 2 - 1 

C 2 (l) = - arccos f -\ + (26) 



= v^ arctan v^-r^- (27) 



2^ + 20(0-1) • (28) 



T/ie Case < /? < 1 
In this case, Fj~ and F~ are given by equation (B12)-(B13) with: 

D 6 (l) = Aarccosh (±) + + ^ + ^) , (29) 

D A (l) = ^arccosh (±) + v 7 !^ (A + _L) , ( 30 ) 

D 2 (l) = iarccosh + (31) 



= J ^arctanJ ^3^- (32) 



(2/3 -Z 2 )*! (/?,/) (/?-/ 2 )7T^ 
* 2 (/M) = ^ + w _ 1} • (33) 

In the particular case of a dominating central mass, the DF is obtained setting (3 = in 
equations (B12)-(B13), eliminating the term containing X 2 , and using equations (29)-(32) with 
lim^ ^i(/5,0 = VT^l 1 . 
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5. Consistency of HH Models 

We can now explore the parameter space of HH models, studying numerically the formulae 
derived in §4. The simplest way to summarize the results is to express the consistency limitations 
in terms of the (normalized) anisotropy radius of the reference component. This is particularly 
indicated because: 1) we know from the introductory discussion (§ 3.2) that the globally isotropic 
component of the HH model is consistent whatever the halo component is (i.e., /; > 0), and 
2) in the case of the OM anisotropy, the anisotropy radius can be isolated in the DF. So, in 
the parameter space (s a ,/x, /3) it is easy to determine the critical value s ac (//, /?) defined as the 
anisotropy radius such that for s a < s ac there exists at least one permitted value of the potential 
for which / < 0. Imposing the positivity of / over all the domain < q < 1, from equation (19) 
one obtains: 



«acGM) = SUp 



Mi) 



(34) 

<?e]o,i[ 



In Fig. 1 different curves s ac are plotted, for a varying halo mass, and for fixed values of f3. 
As expected no model permits a smaller value of the anisotropy radius than that derived from the 
necessary condition. The intermediate solid line is the value of the exact limit on the anisotropy 
radius for the H model, s a ~ 0.202. The upper solid curve is the plot of the WSC in the case of a 
central BH of normalized mass fi. As described in §3.2 this curve approaches asymptotically the 
value l/\/2. 

The dashed and dotted lines represent the lower bound on the anisotropy radius for various 
[3. Clearly all the families of curves for /i — ► converge to the value required by the H model. 
The first result is that the critical anisotropy radius for each model is not strongly dependent 
on the halo mass, and it is always contained between the value given by the necessary condition 
and that obtained from the SSC applied to the H+BH model. The second result is that all 
models with < [3 < 1 (short-dashed lines) have a critical anisotropy radius higher than the H 
model: as expected a model with a concentrated halo cannot sustain also too much anisotropy. 
The long-dashed line represents the case of a central BH: anisotropy radii higher than the values 
represented by this curve can be assumed independently of the halo structure and mass. The 
third result is that all models in which the halo is more extended than the reference component 
can have a slightly smaller anisotropy radius than the H model, but in this case the effect is much 
less stronger than for a (3 < 1 halo: as already found in CP92, the most diffuse component of 
a multi-component system is also the most "delicate" concerning the consistency. Finally, note 
how all the curves become flat asymptotically. Their limiting value cannot be obtained directly 
using equation (19) for very high fi values, due to the numerical precision loss of the intervening 
functions. On the contrary a very accurate analysis can be performed using the asymptotic 
functions given in §4.3, and the result is plotted in Fig. 5. 
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5.1. Stability of HH Models 

As important as the discussion on the consistency of HH models is the location of the region of 
stability for such systems. A complete stability analysis is beyond the task of this work, requiring 
N-body simulations or normal mode analysis, but some interesting conclusions can be equally 
derived, at least for the globally isotropic case. In fact in this case powerful theorems are at our 
disposition: for stability against both radial and nonradial perturbations it is sufficient (but not 
necessary!) that the system DF (in our case the sum of the DF of the halo and of the reference 
component) is an increasing function of the relative binding energy £ (see, e.g., BT87, p. 296-307, 
|Fridman and Polyachenko 1984 , p. 152-163). 



In the performed numerical exploration of the parameter space (p, (5) all the DFs with 
(3 > are monotonically increasing functions of the parameter q, as shown for some particular 
cases in Fig. 2. Changing the parameter q of each component to £ using equation (17), and 
taking the derivative shows that the globally isotropic HH models are stable. Unfortunately the 
sufficient condition above cannot be applied to the globally isotropic case with a central BH: the 
condition df jd£ > is not verified, due to the convergence of the DF to 0, and so only numerical 
investigations can answer this interesting problem. 

For anisotropic systems the situation is more complex, for the lack of general theorems. In any 
case some hint can be obtained by the empirical requirement that £ = 2K ra d/-Ktan^£c = 1.7 ±0.25 



( Fridman and Polyachenko 1984 , p. 235) in order to avoid the radial orbit instability, where i^ ra d 
and K tan are the total radial and tangential kinetic energies of the system, and the exact value of £ c 
is model dependent. A second complication with the previous criterion arises form the fact that a 
generalization to multi-component systems is not obvious, as shown by Stiavelli and Sparke (1991). 
In the same paper the authors show with the aid of N-body simulations that the presence of a halo 
does not change very much the situation with respect to the one-component model. For HH models 
with OM anisotropy the evaluation of i^ ra d an d K ta _ n can be done analytically, but the resulting 
formulae for (3 > are complicated. For these reasons I show here only the trend of the parameter 
£(s a , p) for the particular case of the H+BH system: the simpler case of an H model is obtained for 
p = 0. From the virial theorem K Tax \ + -ftTtan = (\U\ + | W|)/2, where U is the gravitational energy 
of the H model, W is its interaction energy with the BH, and so £ = 2/[(|C/| + |W / |)/2i ; r ra d — !]• 
Choosing as the energy normalization ?7n = M^n, it results U = — 2-7T J p^r 2 dr = —U^/Q, 
W = 47rjr 2 p(d^ h /dr)dr = -U m , and K rad = 2ir J pa 2 r 2 dr = U N (K? ad + pK™), with 

- H l + 8^ + 23 Sa + 12 S 6 ^ a (_l + l 6g 2 + l 5s 4 + 6g 6) 

rad 12(1 + s 2 ) 2 + 24(1 + s 2 r 

, / x ^ \ sJ7 + 8s 2 + 3si) it . . 

In «a> a \ a 2 3 &> + - + arctan( Sa ) 
L 3(l + s£) d 2 

sj[K\n{l + si) + sM-sl2M2)] 

a . ( 35 ) 
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and 



12^ + 1 7T£a(2 + 3 f |) |^a(5 + 6 S 2) 

+ 2(1 + sl) +ln ^ Sa 
,(1 + 6sl)[7rln (1 + si) + s^(-sl 2, 1/2)] 



TT 



(36) 



4 

I have used for of (r) required in the integration the expression given in Ciotti et al. (1996) [eqs. 



(Al)-(A3)], and the function $ is the Lerch transcendent (see, e.g., Erdely et al. 1953] vol.1, p. 27). 



In Fig. 6 the function £ is plotted for various \i values, and the asymptotic flattening to unity 
for increasing isotropy is evident. The first comment is that the stability criterion requires, also 
in the conservative hypothesis of a very high value of £ c (i.e. ~ 2.5), minimum anisotropy radii 
appreciably larger than those obtained form the consistency analysis. For example, the minimum 
anisotropy radius permitted for an H model is ~ 0.93, much higher than the value required by 
the simple consistency, ~ 0.2. So, it is likely that the more radially anisotropic HH models with 
positive DF are with any probability radial orbit instability prone. Another important comment 
concerns the relation between the relative distribution of the halo and of the reference component. 
As can be seen from Fig. 6, haloes more extended than the reference component make the model 
more unstable; on the contrary, more stable systems are obtained for more concentrated haloes. 
This is easily explained looking at the radial trend of the velocity dispersion as a function of the 
anisotropy radius. For extended haloes, the velocity dispersion is mainly increased in the outer 
parts of the model, where also the orbits are strongly radial, and this correspondingly increases 
the £ value; the opposite happens for very concentrated haloes. Particularly explicit are the two 
lines in Fig. 6 referring to a central BH: in this case the velocity dispersion is essentially increased 
only at the model center, where isotropy is nearly realized, and the stability indicator remains low 
also for small values of r a . But in this case probably the global indicator £ looses its meaning, due 
to a very strong decoupling between the central and the outer parts of the model. 



6. Conclusions 

In this paper an extensive analytical investigation of two-component spherical galaxy (or 
cluster) models, made of the sum of two Hernquist density distributions with different physical 
scales, is carried out. A variable amount of orbital anisotropy is also allowed in both components. 
These models, characterized by a power-law density profile in their central regions - both in the 
visible and in the dark matter distribution - reproduce the main properties of early-type galaxies 
as revealed by Hubble Space Telescope observations, and also of the dark matter distribution 
as obtained in recent N-body simulations. The main results presented in this paper can be 
summarized as follows: 

1 . The analytical expression for the DF of HH models with general OM anisotropy is presented 
and discussed, even for the particular case of an Hernquist model with a central BH. The 
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special case of a dominant dark halo is also discussed, and it is shown that under this 
assumption the DF can be asymptotically expressed using just elementary functions. In the 
case of global isotropy the analytical expression for the differential energy distribution of 
both components is obtained. Some velocity sections of the DF are shown and discussed. 

2. The necessary and sufficient conditions that the model parameters must satisfy in order 
to correspond to a consistent system (i.e., a system for which each physically distinct 
component has a positive DF) are analytically derived using the method introduced in 
CP92. It is proved that globally isotropic HH models are consistent for any mass ratio and 
core radii ratio, even in the special case in which the "halo" reduces to a BH. In the case of 
a central BH and of variable anisotropy for the host system, the analytical expression for a 
minimum anisotropy radius as a function of the BH mass is given. 

3. The region in the parameter space in which HH models are consistent is subsequently 
explored using the DF. The main result is that the presence of a massive halo does not 
affect significantly the maximum anisotropy that can be sustained by a consistent model. 
It is shown that the presence of a halo with a core radius larger than that of the reference 
component allows a slightly higher degree of anisotropy with respect to the one-component 
Hernquist model. On the contrary, a halo with a smaller core radius imposes a larger value 
for the minimum anisotropy radius than that proper of the H model. The most restrictive 
case is that of a central BH. In any case, for a given core radius of the halo there is a lower 
limit to the minimum anisotropy radius that approaches an asymptotic value for a dominant 
halo mass. 

4. Finally, it is proved that isotropic HH models are stable, except for the case of a central BH, 
when no conclusions can be drawn. For anisotropic models the stability parameter against 
radial orbit instability is briefly discussed, and it is shown that with high probability the 
most anisotropic consistent HH models are unstable. 



I would like to thank Giuseppe Bertin and Silvia Pellegrini for helpful comments and 
discussions. This work has been partially supported by the Italian Ministry of Research (MURST). 



A. Consistency Requirements 

In this Appendix I apply the CP92 method to the HH models, as discussed in § 3. The 
necessary condition for the H model imposes through equation (4) a limitation on the anisotropy 
radius: 

(3tt-2)(l-*)3 - 
a " (4-3^)^ 2 ' - - v ; 
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The requirement is that si is larger than or equal to the maximum of the function on the r.h.s.. 
This is reached at 

*m = \ - ^§ c 0.77, (A2) 
that, after back substitution in equation (Al), gives equation (13). 
The SSC applied to the H model gives the following inequality: 

2 3(5* - 2)(1 - *) 3 ~ /t , 

s 2 a > — i= ^ < * < 1. (A3) 

a " (15* 2 -39^ + 28)^ 2 _ _ V ; 

After differentiation, discarding the two complex conjugates roots of the resulting cubic equation, 
and defining 

/ 1336 V35887965 \ 1/3 , AJ , 
*° = [m& + 16200 J ' (A4) 

the maximum of the r.h.s. of equation (A3) is reached at 

23 217 

*-= 18 + 1620S-^ - 52 - <A5) 
that after back substitution in equation (A3) gives equation (14). 

The WSC applied to the H+BH model is obtained from equation (6) with b = — 1. The 
inequality to be verified for a given \x is: 

> _3*4 _ 10*3 + 6(2 + ^2 _ 6(1 +Jl) * + {1 + fl) o< _ <l 
(3* 2 - 8* + 6 + 6^)* 2 - - 

Taking the derivative of the r.h.s. we are left with the discussion of a rational function whose 
denominator is strictly positive, and the numerator is a polynomial of fifth degree, that fortunately 
can be factorized in a term of second degree strictly positive, and the cubic: 

c M (*) = -* 3 + 4* 2 - 6(1 + ii)if + 2(1 + fj). (A7) 

Observing that c^(0) = 2(1 + jj) > and c^(l) = —(1 + 4fi), there is at least one solution of 
c M (*) = between and 1. This solution is also the only one, because the determinant of the 
cubic c M 

_ 4(l + / u)(216/x2 + 99^+ll) 
ii 27 ' 

is positive V/i > 0: from the theory of algebraic equations, c^(*) = admits one and only one 
real solution. Finally, considering the sign of we proved that the solution *m corresponds to a 
maximum. Defining 



V3(l + //)(216// 2 + 99/i + llJ _17 
9 ^27 



1/3 

(A9) 
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one obtains 

*mM = |-^^ + ^, (MO) 

that after substitution in equation (A6) gives equation (15). When /j, = 0, the value at which the 
maximum is reached can be calculated directly from equations (A9)-(A10), obtaining: 



4>o 

and 



V33 _ 17 
~~9 27 



1/3 

(All) 



= \ ~ ^7- + - 0.456, (A12) 

3 y<po 

that after substitution in eq. (A6) with \i = gives eq. (16). 

The application of the WSC to the globally isotropic HH models is more complicated. After 
having computed the derivatives in equation (6), we have to investigate the positivity of a rational 
expression, whose denominator is strictly positive Mb > —1, b / 0, V/z > and < VP < 1. The 
numerator of this function factorizes in a strictly positive function and in the polynomial: 

A«%(tf) + (3^ 2 - 8* + 6)(1 + b^f > 0, (A13) 

where 

c b (*) = 66^ 3 + 3(1 - 56)§ 2 + 2(56 - 4)* + 6. (A14) 

It is trivial to show that the second addend of equation (A13) is strictly positive. I prove 
now that also c& is positive for b > — 1 and for < VP < 1: so the WSC is satisfied for any 
choice of (/u;/?) and the globally isotropic HH models are consistent. First of all, Q,(0) = 6 and 
0,(1) = 6+ l = /3>0, and so after excluding the presence of roots of q, in the interval [0,1], 
the proof is obtained. This can be accomplished using the classical Sturm method, i.e., counting 
the differences in the number of variations between and 1 in the Sturm sequence Sb associated 
to Cb(*S>) and discussing it as a function of the parameter b. But a faster proof can be obtained 
as follows. Since q, does not change sign between and 1, only an even number of roots can be 
contained in the interval. If we show that equation (A14) admits only one real solution then the 
proof is obtained. The discriminant of q, 

Q3& + Mt + «)(> + If 

29166 4 V ; 

is strictly positive, and so the cubic equation admits one and only one real solution, that is 
necessarily placed outside the interval [0,1], and q,^) > for < \P < 1. 
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B. The DF for HH Models 

In this Appendix the main analytical steps required for the determination of the DF are 
described. As discussed in §4, the first step is to change the integration variable from the total 
potential to the potential of the reference component. After a normalization, equation (2) becomes: 



HQ) 



/ N dF[q{Q)\ 



V8ir 2 dQ ' 

where the relation between q and Q is given by equation (17), and 

rv dg d^f 



F(q) = 

From equations (11) and (17) 

§t(q) - *t(*) 



M Jy T (q) - * T (*) 



(Bl) 



(B2) 



(g- + l + n 
— , o = -k 



1 + 6* 



(B3) 



It is shown in the next subsections that F is actually a function of q only through I defined by 
equation (18). Unfortunately its form changes depending on the sign of b = f3 — 1, and it is 
convenient to separate the discussion of the cases (3 > 1 and < (3 < 1. 



B.l. The Case (3 > 1 

This is the simpler case, for which b > and 1 < I 2 < (3. After the change of variable 



equation (19) becomes: 



F+(l) 



2 r l 



*■<«) + 5# 



dt 



6 5 /2 J x 

After a conversion in simple fractions it results: 

K a {t) = -3t 6 + 2(2 + (3)t 4 - (1 + 2(3)t 2 , 

and 



Ki(t) = -3i 6 - 2((3 - 4)i 4 - ((3 2 -4(3 + 6)t 2 



(0-1)4 (3(p-iy- 



+ 



(B4) 



(B5) 



(B6) 



(B7) 



(3-t 2 (/?-t 2 ) 2 ' 

Using the nomenclature of Byrd and Friedman ( Byrd and Friedman 1971] , hereafter BF71) we 



have: 



2g_ 

b 5/2 



fin o/« ,M4- ,„2 (/3-l) 4 T , , /3(/3-l) 4 



-3Z b C 6 - 2(/3 - 4)PC 4 - (/? 2 - 4/3 + 6)/ 2 C 2 



/3-Z 2 



(/3-Z 2 ) 2 



Vo 



(B8) 
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and for the anisotropic part, 



2g_ 

65/2 



[-3l 6 C 6 + 2((3 + 2)/ 4 C 4 - (1 + 2(3)l 2 C 2 ], 



where the functions C m (4>,k) and V m (4>, a 2 , k) are given in Appendix D, and 



k 



V 



S 2 + l 



, <fi = arccos 



a 2 = 



i 2 -/?' 



(BF71, p.48). 



(B9) 



(BIO) 



B.2. The Case < (3 < 1 

In this case f3 < I 2 < 1. After the change of variable given in equation (B4), equation (19) 
becomes: 

dt 



F-(l) 



2 



(Bll) 



For the isotropic part: 



|6| 5 / 2 

and for the anisotropic part 



3l 6 D 6 - 2(13 - A)l 4 D 4 - ((3 2 - 4/3 + 6)l 2 D 2 



x x + 



W-i 2 )'- 



-Xo 



|6| 5 / 2 ' 



-3/ 6 D 6 + 2(/3 + 2)/ 4 D 4 - (1 + 20)l 2 D 2 ], 



where the functions D m ((j), k) and X m (cj), a 2 , k) are given in Appendix D, and 



1 



8 2 + l 2 ' 



arccos (I), a 2 = 



f3-l 2 ' 



(B12) 
(B13) 

(B14) 



(BF71, p.45). 



The particular case of a central BH can be derived from the previous formulae. In this case 
b = — 1 and so < I 2 < 1. The coefficients of the special functions are obtained from equations 
(B12)-(B13) for (3 = 0, and their arguments are still given by equation (B14). The only problem 
is in the isotropic part: the coefficient of X 2 is zero and lim Q ,2^ X\(cj), a 2 , k) = C 2 ((f>,k). The 
anisotropic part is obtained by direct substitution of (3 = in equation (B13). 



C. The Density of States 

The formula (4.157b) of BT87, after normalization of radii and energy to the core radius and 
central potential of the reference component can be written as: 

_ _ rsM(£) „ rz z 

G(£) = J s 2 ^ T (s)-£ds, (CI) 
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with 'I't('Sm) = £■ Changing the variable of integration from the radius s to the relative potential 
of the reference component, the result is: 



"(1 - § 



i2 



G{q) = J V*tW " *T(g) ^ 4 d*, (C2) 

where the relation between g and £ is given by equation (17). Again as for the DF it is convenient 
to discuss separately the two cases (3 > 1 and < (3 < 1. 



C.l. The Case > 1 



Changing the variable of integration as in equation (B4), from equation (C2) 

'VP G(t)dt 

y/(P+P)(!?-F) 

where 

2/ 3 _ 4 + Z 2 - J 2 (^_/2 )(3 _ 2/ 3 )+/?(/ 3_ 6)+6 _ /U 

= 1 + r~fi + {T-^y + 



(1 -/?)[(/? -3)(<5 2 -Z 2 ) + 2(/3 + ^- 2)] (1-^(1-/3)2(1 + ^) 

(1 - i 2 ) 3 (1-i 2 ) 4 ' 1 J 



After integration, 



G+(0 _ 2/3 - 4 + Z 2 - J 2 ^ _ f 2 )( 3_ 2/ 3) +/?(/?- 6 ) + 6- 

2^ " l3]2 X l+ (7372)2 ^ + 

(1 -/?)[(/? -3)(<5 2 -Z 2 ) + 2(/3 + M - 2)] y (1-OTl+g) 

^372)3 ^ + ^2)3 *4, (C5) 

where the functions X m are given in Appendix D, and their arguments are 



52 A ( 1 \ 
. if) = a.rrrns 



1 



9= 7FTP' k= 5^TP> ^ = arccos^— j, a = (C6) 



C.2. The Case < (3 < 1 

In this case, 



where Q{t) is given in eq. (C4). After integration, 

G~(l) 2 p-A + P-5\ r (a2-;2 )(3 _ 2/ 3 )+/?(/? _ 6) + 6 

= ^0 H ; P5 ^ 1 H 7^ 7775 ^2 + 



-2|6| 1 /2 5 " u i_/2 ^ (1-/2)2 

(1 -/?)[(/? -3)(J 2 -Z 2 ) + 2(/3 + j u- 2)] (1-/3)2(1 + ^ 

^2)3 V * + (1 _ /2) 3 ^, (C8) 
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where the functions V m are given in Appendix D, and their arguments are 

The case of the central BH is obtained from the previous two equations with the substitution 
(3 = 0. 



D. Elliptic and Jacobian Functions 

Here I shortly summarize the notation for the special functions used, and their mutual 
relations, following BF71. The zero-th order functions are: 

C (<f>, k) = D (<f>, k) = V (<f>, k) = F(<f>, k). (Dl) 

The m-th order functions are given by: 

C^k) = mk) -£- k2)C \ (D2) 

n , \ 2m(2fc 2 - l)C 2m + (2m - 1)(1 - fc 2 )C 2m „ 2 + sn(^)dn(n)cn 2m -» m - 
C 2m+2 (</>, k) = {2m + l)k 2 ' 

(BF71, eq. 213.06-312.05). 

V^a^k) =U( ( f ) ,a 2 ,k), (D4) 

, 2 , , a 2 E(<f>, k) + (k 2 - a 2 )V + (2a 2 k 2 + 2a 2 - a 4 - 3A; 2 )Vi 
1/2(0, a; ,k) = 



2(a 2 - l){k 2 -a 2 ) 
a 4 sn(ii)cn(-u)dn(u) 
2(a 2 -l)(fc 2 -a 2 )[l-a 2 sn 2 (u)]' 



w 1 ± 2 t \ (2m + l)fc 2 y m + 2(m + l)(a 2 fe 2 + a 2 -3fe 2 )F m+ i 
F m+3 (<^,fc) - 2(m + 2)(l-a2)( ib 2_ Q 2 ) + 

(2m + 3)(a 4 - 2a 2 A; 2 - 2a 2 + Zk 2 )^^ 



2(m + 2)(l-a 2 )(fc 2 - a 2 ) 
a 4 sn(tt)cn(n)dn('u) 



2(m + 2)(1 - a 2 )(fc 2 - a 2 )[l - a 2 sn 2 (u)] m + 2 ' 

(BF71, eq. 213.11-336.03). 



(D5) 



(D6) 



D 2 ^ k) = (l-m-^fe)+tn(n)dn ( ,) ; (D?) 

1 rC 

M 2m(l - 2fc 2 )Z? 2m + (2m - l)^^^ + tn(n)dn(n)nc 2m (u) 

^ 2m+2 ^' fc) = (2m+l)(l-fc2) ' (D8) 
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(BF71, eq. 211.09-338.05). 

1 m l 
X^a^k) = ]>> 2 - iy m - Vj , (D9) 

(BF71, eq. 211.14-338.04). 

In the previous expressions, the elliptic integrals of first, second and third kind are expressed 
as functions of the phase (ft and the modulus k: 

F(<P,k)= / : (D10) 

Jo VI - k 2 sin 2 

E{<j>,k) = Jl-k 2 sin 2 6d6, (Dll) 
•/ o 

U(<P,a 2 ,k)= / ; (D12) 

Jo (I -a 2 sin 2 Q)\J\-k 2 sin 2 6 

The relation of the elliptic integrals with the Jacobian functions for given (4>, k) are: 



u = F((fi,k); cn(u) = cos(0); sn(u) = sin(0); dn(n) = y 1 — A; 2 sn 2 (u), (D13) 

and 

tn(n) = ^; nc(u) = (D14) 
cn(n) cn(it) 
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Fig. 1. — The limits on the anisotropy radius for the consistency of the reference component, for 
various (//, f3) values. The lowest solid line is the necessary condition. The middle solid line is 
the minimum anisotropy radius for the H model. The dotted lines are the limits in the presence 
of a more extended halo, with (3 = 2,5, 20. The short-dashed lines correspond to a halo with 
(3 = 0.5, 0.1, and the long-dashed lines to the case of the central BH. The dot-dashed and solid line 
are the SSC and WSC respectively, for the case of the central BH. 

Fig. 2. — Upper Panel: the DFs for the reference component of the HH models in case of global 
isotropy. The solid line is the DF of the H model. The dotted lines correspond to the DF in 
presence of a halo more extended than the reference component, with (3 = 5, and \i = 5, 10, 20. The 
short-dashed lines refer to a more concentrated halo, with (3 = 0.2, and with masses \x = 0.1, 1, 10. 
Finally, the long-dashed lines correspond to a central BH, of masses /i = 0.01,0.1, 1. Lower Panel: 
the DFs for the reference component of the HH model in the case of radial anisotropy with s a = 1. 
The various curves correspond to the same parameters as in the isotropic case. 

Fig. 3. — Radial (solid lines) and tangential (dotted line) normalized velocity sections for the H 
model, with s a = 1. The numbers near the lines give the normalized radius at which the section is 
shown. The dashed line is the section at r/r c ~ 0.01 in the case of a central BH with /i = 0.1. 

Fig. 4. — The differential energy distribution for the reference component in case of global isotropy. 
The solid line is the dM/d£ of the H model. The other lines correspond to the same parameters 
described in Fig. 2. 

Fig. 5. — The minimum value for the anisotropy radius in case of a dominant halo. Note the 
qualitative difference moving from a less extended halo to a more extended halo than the reference 
component. 

Fig. 6. — The value of the stability parameter for the reference component as function of the 
anisotropy radius. The solid line refers to the H model, and the dotted lines to H+BH models: 
the numbers are the ^'s values. The other lines are the value of the stability parameter for various 
values of (/j,; (3): circles (1,0.1), triangles (10,1.5), and squares (10; 5). 
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